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We simulate the growth of a dual species Bose-Einstein condensate using a Gross-Pitaevskii equa- 
tion with an additional gain term giving rise to the growth. Such growth occurs during simultaneous 
evaporative cooling of a mixture of two gases. The ground state of a dual condensate is normally 
either a miscible mixture, or an immiscible phase with two spatially separated components. In a 
cigar trap the ground state typically consists of one component in the center, and the other compo- 
nent flanking it. Our simulations show that when the condensates are formed in a cigar trap and 
the mixture is phase separated, then the final state upon the end of the growth is generally far from 
the true ground state of the system. Instead it consists of multiple, interleaved bubbles of the two 
species. Such a pattern was observed recently in an experiment by Wieman's group at JILA [l|], 
and our simulations are in good qualitative agreement with the experiment. We explain the pattern 
formation as due to the onset of modulation instability during growth, and study the dependence 
of the final state pattern on various parameters of the system. 



I. INTRODUCTION 

Like common fluid mixtures, dual species Bose- 
Einstein condensates (BECs) may exhibit both miscible 
and immiscible phases. In the miscible phase, the two 
species form a homogeneous solution, while in the im- 
miscible phase, they separate spatially 0,3- I^ ^ BEC, 
the occurrence of these phases is controlled by the ratio 
of the inter-component interaction to that of the intra- 
component interactions. 

The ground state structure of the immiscible phase is 
typically that of a ball and shell where one species forms 
a shell around the other, or that of two condensates lying 
side by side. In a quasi-lD trap, the structure can be of 
one component flanked by the other. The exact structure 
depends on the shape of the trap, the size of each com- 
ponent, and the interactions between them. In general, 
the two components tend to minimize the surface area 
that separates them, as this surface contributes surface 
tension energy [1 i, i, S • 

By starting from the ground state of a miscible phase 
and increasing the inter-component repulsion (relative to 
the intra-component interactions), it is possible to make 
a transition from the miscible to the immiscible phase. 
Such a change in the interactions can bring about a mod- 
ulation instability where the two components separate 
locally and form a pattern of interleaved bubbles. Alter- 
natively, such a pattern emerges when half the population 
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of one pure component is suddenly transfered to another 
one, and the two components separate [a [3, [lOf| . 

A recent experiment p^] on a mixture of ®^Rb and 
^^Rb in an elongated cigar trap, observed an intriguing 
modulated pattern upon the end of evaporative cooling 
(shown in figure 4 of that reference). The pattern is that 
of separated and interleaved "bubbles" of the two com- 
ponents, with up to 3-4 bubbles from each component. 
This behavior is seen when the inter-component repulsion 
is strong compared to the intra-component interactions, 
i.e, when an immiscible phase is expected.. These spe- 
cific patterns were observed even under conditions when 
the interactions strengths were kept constant during and 
after the evaporative cooling. Nor is there population 
transfer involved, as the two components are different 
isotopes. Following evaporative cooling one may then 
expect to obtain a binary condensate in its ground state. 
But the observed modulated pattern does not look like 
the ground state of the system. 

In this paper we suggest that the emergence of this pat- 
tern has to do with the dynamics of the formation of the 
binary condensate during the evaporative cooling. When 
the condensates are formed they are initially very small. 
Their mutual interaction (proportional to the densities) 
is then weak and we expect them to coexist spatially, 
even when the predicted phase for a Thomas- Fermi (large 
number of atoms) regime is immiscible. As the conden- 
sates grow, spatial separation is triggered by the growing 
inter-component repulsion. One may expect this to bring 
about a modulation instability giving rise to the final ob- 
served spatial pattern. 

Our aim in this paper is to demonstrate this mecha- 
nism, predict the conditions for the formation of the spa- 



tial structure, and determine the factors that influence 
the number of bubbles formed. The paper is organized 
as follows: in section (jll| we give the theoretical back- 
ground: we present the Gross-Pitaevskii gain model for 
a double condensate, and review the theory of modula- 
tion instability mechanism. In section (IIII|) we present 
and discuss the results of our numerical simulations. We 
summarize our conclusions in section (IIVI). 



II. 



BINARY CONDENSATE GROWTH 
MODELING 

A. Gross-Pitaevskii Gain model 



We are interested in the dynamics of the formation and 
growth of the binary condensate during the evaporative 
cooling, in order to understand the final spatial pattern 
that emerges. The description of BEG growth is a chal- 
lenging theoretical problem. In principle it is possible to 
perform first-principles numerical simulations of the rele- 
vant equations [ll|, Il4|- However, a great deal of physical 
insight can be obtained from a simpler model. In the first 
instance, the thermal cloud is much more dilute than the 
condensate. Therefore, it is plausible to treat the sys- 
tem as a two fluid model, where the thermal component 
is assumed to be an ideal gas and its interactions with 
the condensate are neglected. If one is only interested in 
the condensate, then the role of the thermal component 
is simply that of a particle reservoir feeding the growing 
condensate. One thus arrives at the Gross-Pitaevskii gain 
model [131 . This is a Gross-Pitaevskii equation modified 
by a linear gain term. For two components, the Gross- 
Pitaevskii gain equations (GPGE) read: 
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Here, ^1 and ^1^2 are the wavefunctions corresponding to 
the two condensates. The normalization of each wave- 
function is taken independently as /dr|^i(r)p = iV^. 
mi and TO2 are the corresponding atomic masses. The 
condensates are assumed to be trapped in axis-symmetric 
harmonic potentials. 
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KS(^, z) = -m, yy + ui^,){z - Z,Y\ , z = 1, 2, (2) 

where u!r,i^z(i) &re the radial and transverse trapping 
frequencies, and Zi the axial center of the traps. We allow 



for the possibility of different axial trapping frequencies 
for the two condensates, as well as for a relative axial 
shift Z2 — Z\ between the centers of the traps, as occurs 
in the experiment [ll. 

The intra-component coupling constants gi — 
Anh ai/rrii are characterized by the scattering lengths 
fli and a2 between atoms of the same species, while 
the inter-component coupling, 512 = 27r?i 012/77112, (with 
— !— = -^ -I- — ), is determined by the scattering length 
ai2, where an atom of component 1 scatters from an 
atom of component 2. The P^'s are the gain terms. 
They contribute to the growth of component i within 
timescale l/Pi. Intuitively, the growth terms arise from 
boson statistics, where the probability to add an addi- 
tional atom to condensate i is proportional to the num- 
ber of atoms already condensed. In practice these are 
phenomenological parameters and may be estimated by 
fitting to an experiment. 



B. Modulation Instability Mechanism 

Following the analysis of Refs. [H S Q we define the 
parameter 



A = 9192/ gl2 - 1 



(3) 



that depends upon the ratio of the single species and 
interspecies interactions. When the number of atoms is 
large enough the system may be treated by the Thomas- 
Fermi approximation, which neglects the kinetic energy. 
Assuming 912 > 0, there are two regimes: A > 0, where 
the two condensates are miscible, and A < 0, where they 
are immiscible due to the interspecies repulsion. 

It is important to bear in mind that this classification 
is correct only when the Thomas-Fermi approximation is 
valid. It does not hold for small, trapped condensates, 
where the kinetic energy terms are significant. In the 
numerical simulations presented in the next section we 
shall show that in this case the condensates can initially 
co-exist spatially even when A < 0. As they grow, spa- 
tial separation sets in. To have a better understanding 
of the dynamics of spatial separation, it is useful to con- 
sider an idealized system, of two homogeneous conden- 
sates in a quasi one dimensional geometry (relevant to 
the experiment [l|), with respective densities n\ and 712, 
co-existing spatially at time i = 0, as discussed in de- 
tail in Refs. [3, [lO|- If A > 0, the system is stable. If 
A < 0, the initial state is unstable, and an initially homo- 
geneous system will evolve dynamically. In this so called 
modulation instability mechanism, the two components 
separate, and the separation occurs in a wave pattern 
with a typical length scale which depends on the rela- 
tive interactions and densities. Regions of over-density 
of one component appear in regions of under-density of 
the other component, and vice versa. Nonlinearity ac- 
celerates the process and brings about complete spatial 
separation. A pattern of interleaved bubbles of the two 
components emerges. The final spatial pattern can be 



very stable, since there is no room for movement to "the 
sides" , and the energy barrier for tunnehng between the 
separated bubbles is very big. 

The analysis of the modulation instability in a dual 
component BEC in a quasi ID geometry was done in Ref. 
[l(| . There, the assumption of equal masses mi — m2 = 
m has also been made, which considerably simplifies the 
analysis. The behavior of the system at small i > can 
be deduced from the excitation spectrum. The eigen- 
modes of the system are characterized by a longitudinal 
wave number k and frequency fi. For each eigenmode, 
the wave number k is shared by the two components [l6| . 
The instability manifests itself in a continuous range of 
modes which have imaginary frequencies. The most un- 
stable mode is that with imaginary frequency of maxi- 
mum absolute value. This is the mode that will grow 
fastest. Its wavenumber kmax is found to be: 
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where b is the radial harmonic oscillator length b = 

— ^, and we assume here ■m^ = mo = to. 

The imaginary frequency of this mode is found to be: 
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FIG. 1: Calculated ground state of a dual (*^Rb/ ^''Rb) con- 
densate in an immiscible phase, in a cigar shaped trap. Shown 
are the axial density profiles. The interaction parameters are: 
as5 — 200 bohr, agr — 99 bohr, ag5_87 — 214 bohr. The trap 
parameters are: radial frequency 2tt x 130 Hz for both com- 
ponents; axial frequency 27r x 2.9 Hz for *^Rb and 2-k x 2.6 Hz 
for *^Rb. The number of atoms is 5 x 10* in each component. 
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and the time scale for growth is Tmax — 27r/G„ 



III. NUMERICAL SIMULATIONS 
A. Mechanism of dynamical phase separation 

We have performed numerical simulations of the 
growth of the dual *^Rb/^^Rb condensate based on 
Eq. ([T]). In correspondence with the experiment ll|, 
all our simulations have the following fixed parame- 
ters. asT — 99 bohr is the scattering length of ®'^Rb. 
^85-87 = 214 bohr is the inter-species scattering length. 
The trap parameters are: radial frequency 27rx 130 Hz for 
both components; axial frequency 27r x 2.9 Hz for *^Rb 
and 27r x 2.6 Hz for ^^Rb. 

We first demonstrate that the true ground state of the 
system in an immiscible phase is indeed quite simple, as 
described above. We computed the ground state with 
a85 = 200, and the other scattering lengths fixed as 
above, so that A < 0. The ground state for this sys- 
tem with 50,000 atoms in each component is shown in 
Fig. 11]). The calculation is 3D, and the figure shows 
the density profile p{z) along the long z-axis, where 
p{z) = J p(x, y, z) dx dy. We have also performed addi- 
tional computations which modeled the effects of a grav- 
itational sag that breaks the cylindrical symmetry and 



slightly shifts the centers of masses of the the two com- 
ponents. In all cases, the ground state was still similar 
in shape to that of Fig. ([1]). 

We now investigate the consequences of a dynamical 
growth of a dual condensate during evaporative cool- 
ing. The exact growth time of the condensates in the 
experiment [l| is not known, but has been roughly es- 
timated to be between 100-500 ms [17]. In our simu- 
lations we begin with the ground state consisting of 10 
atoms of each component, and let the condensates grow 
for 300ms. We choose gain terms Fi = r2 = 14.2sec~^, 
giving rise to a final population of 51,000 atoms in 
each component. The number of atoms Ni{t) of com- 
ponent i at time t can be calculated from: Ni(t) = 
iVi(0)e^'"'*. We assume that after 300ms full condensa- 
tion has been achieved and the growth is terminated. 
The Fi's are then set to zero and further dynamics 
are atom number conserving. Fig. ^ shows stages 
in the growth of the dual condensate under these con- 
ditions. Animations of this evolution are available at 



http://grizzly.colorado.edu~bohn/movies/bubbles.htm 



In the initial state (panel i = 0), the number of atoms 
is very small so that the mean field interactions are weak. 
Each condensate's shape is essentially Gaussian and they 
are almost perfectly overlapping. As the density grows 
(due to the imaginary gain term), the interactions kick 
in (around 130ms, compare panels t — 120 and t — 180), 
leading to phase separation at the center, where density 
is largest. After the end of the growth (at 300ms, panel 
t — 300) the modulation pattern is observed to be quite 
stable, apart from a breathing motion, with the compo- 
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FIG. 2; Dynamical growth of a dual BEG with the same interaction and trap parameters as in Fig. ([T]). The initial state is the 
ground state with 10 atoms in each component. 
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FIG. 3: The final spatial pattern resulting from the dynamical 
growth of the dual condensate as simulated in Fig. (Q . 



nents almost fully separated and forming multiple inter- 
leaved bubbles, as shown in Fig. ^. Note the increasing 
vertical scales, which account for growing particle num- 
bers. 

The onset of the immiscibility is triggered by the 
emergence of the modulation instability. The natural 
wavelength of the modulation instability corresponds to 



the unstable mode with maximum imaginary amplitude. 
This is 2Tr /kmax where k^ax is given by Eq. ([5]); for the 
purpose of using this equation, which assumes species 
with equal masses, we use the average mass of the two 
isotopes of Rb, neglecting their small mass difference. 
Also, for our non-homogeneous system, we use in it the 
maximum density of each component. As the density of 
the condensates grows, their spatial extent increases and 
at the same time the modulation instability wavelength 
decreases. The decrease of the modulation instability 
wavelength with time is shown in Fig. ([5]) . Around 60ms 
after the beginning of the growth, the modulation insta- 
bility wavelength becomes comparable to the size of the 
condensate as measured by the r.m.s of the density pro- 
file of the combined condensates. The onset of dynamical 
separation is seen somewhat later at 130ms. This lag in 
time may be due to two factors. Firstly, the MI wave- 
length Eq. ^ assumes a homogeneous condensate. In 
other words, it applies when the local density approxi- 
mation (LDA) is valid. On the other hand, the number 
of atoms in each component at t = 130ms is only about 
400. Thus, the LDA approximation is not yet reliable. 
Secondly, as we discuss below, there is a time scale as- 
sociated with the growth of the unstable mode, and this 
leads to a lag between the initial seeding of the instability 
and its growth to an observable size. Nevertheless, the 



crossing of the curves of the r.m.s size of the condensate 
and the modulation instabihty wavelength indicates, to 
within a factor of two in time, the onset of the dynamical 
instability. 

As condensation proceeds and the density continues to 
grow, the spatial extent of the condensates increases and 
at the same time the modulation instability wavelength 
decreases further. This leads to the continuous sepa- 
ration of the two components into smaller and smaller 
"bubbles", as can be seen in Fig. ^. The final size of 
individual bubbles is determined by the point at which 
essentially complete spatial separation between the two 
components is achieved, so that no more sub-divisions 
are possible. The exact point at which this happens de- 
pends strongly on the non-linear dynamics and is thus 
difficult to estimate solely from the linear small pertur- 
bations analysis which is behind Eq. ([5]). Typically, the 
final modulation instability wavelength is much smaller 
than the size of the observed bubbles, but this has no 
effect since at this point the bubbles are already almost 
purely separated components. 

Regarding the temporal evolution, we can compare two 
time scales. One is the time scale over which the conden- 
sate's population grows, and the other is the time scale 
for the growth of the most unstable mode. This modu- 
lation instability time scale depends on the density and 
thus is itself time-dependent. Fig. (O shows the evolution 
of the modulation instability time scale Tmax = q ^ 
as function of the simulation time, where Gmax is given 
by Eq. ([5]). As the density of the condensate grows, the 
timescale for the growth of the most unstable mode is de- 
creasing exponentially. Around 130ms the time scale of 
growth of the unstable mode becomes comparable to the 
time scale of growth of the condensates themselves. As 
we have seen, this coincides with the point where dynam- 
ical phase separation begins to be observable in the simu- 
lation. With increasing simulation time, the modulation 
instability growth time scale becomes very short, leading 
to the rapid appearance of bubbles between 180ms and 
240ms in Fig. ^. 



B. Factors affecting the spatial pattern formation 
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FIG. 4: The modulation instability wavelength (solid line) 
corresponding to the simulation of Fig .([2J, as computed from 
Eq. ((5]|. This is compared with the r.m.s size of the total 
density of both condensates along the trap long axis (dashed 
line). 
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FIG. 5: The timescale for growth of modulation instability 
corresponding to the simulation of Fig. ([2| , as computed from 
Eq. ([5]). The horizontal axis is the simulation time, while 
the vertical axis is the time scale for the growth of the most 
unstable mode. This timescale is Tmax = 1-e jGmax- 



We explored some of the factors that affect the final 
spatial pattern upon the end of evaporative cooling. As 
noted above, the condition for phase separation in the 
Thomas-Fermi limit is A < where A is defined in 
Eq. ([3]). We find that having a negative A with larger 
absolute magnitude leads to the dynamical formation of 
fewer bubbles. For example, Fig. ([6]) shows the modula- 
tion pattern after lOOOms with the only difference from 
Fig. ^ that 085 has been reduced from 200 to 81 bohr, 
so that A, which for the previous case equals —0.57, is 
now —0.82. The number of separated bubbles is clearly 
reduced, from 25 in the first case to 17 in the second. 
The physical reason for this effect is as follows. When 
A is more negative (larger absolute value), the dynami- 



cal phase separation begins earlier in the growth process 
when the size of the condensates is smaller. Since the size 
of the condensates is smaller, they divide into smaller 
number of bubbles. Saturation (that is, complete spa- 
tial separation of the two components) is also achieved 
earlier, and no further sub-division occurs at later times. 
On the other hand, with smaller absolute value of A the 
condensates have a chance to increase there size before 
the onset of phase separation. Thus, eventually, when 
phase separation takes place, a larger number of bubbles 
can be formed. 

An even more dramatic effect is seen when the centers 
of the traps of the two components are shifted one from 
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FIG. 6: The final spatial pattern resulting from the dynamical 
growth of the dual condensate, with scattering parameters 
ags ~ 81 bohr, as? = 99 bohr and ags-g? = 214 bohr; to be 
compared with Fig. ([HJ. All other parameters are the same 
as in Fig. ([3]). 



the other along their long axis by a small distance, as 
may occur due to residual gravitational sag. This can 
occur when the trap is slightly tilted off the horizontal 
to the earth, as apparently occurred in the experiment 
[l|. Even if the tilt is slight, the effect can be large 
in the axial direction due to the low and different trap- 
ping frequencies for the two components in this direction. 
This leads to asymmetry in the spatial pattern formation. 
Fig. ([7]) shows an asymmetrical pattern obtained with a 
dual condensate with the same interaction parameters 
as of Fig. ([3]), but with the centers of the traps of the 
two components shifted 1.7 /im from each other. Fig. ([8]) 
shows the time dynamics leading to this pattern. The ini- 
tial shift is clearly seen in the partial overlap between the 
"seed" Gaussian of the two components. As the density 
grows, the inter-component repulsion causes the centers 
of mass of the two components to move farther apart. At 
the same time, modulation instability and bubble forma- 
tion occurs where the two components overlap. Due to 
the reduced spatial overlap (compared to the case with 
no shift of the traps' centers), the final number of bubbles 
is much reduced. 

We also note that the modulation instability mecha- 
nism which is at the heart of the spatial pattern forma- 
tion requires a quasi one dimensional geometry. In a fully 
3D trap (with aspect ratio close to one), even if "bubbles" 
are seeded, they would quickly merge and coalesce. This 
is because a spatial pattern which contains many bub- 
bles would have higher energy then one with continuous, 
bulk configuration with no bubbles, due to the contribu- 
tion of surface tension energy [3]. But, in the quasi ID 
trap there is a high energy barrier for the merging of two 
bubbles of one component separated by a bubble of an- 
other component. This is essentially a tunneling process 
which is highly suppressed. Our numerical simulations 
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FIG. 7: The final spatial pattern resulting from the dynamical 
growth of the dual condensate, with parameters identical to 
that of Fig. ([3|, except that there is a shift of 1.7 /im between 
the centers of the traps of the two components. 



support this picture. For example, we find that in a trap 
with aspect ratio A = ujr/ojz = 10 (compared with A = 45 
in the simulations discussed above) it is quite hard to see 
bubbles, and at most four to five are formed. 

Another factor in the growth dynamics concerns the 
start of condensation of the two components. When the 
atom numbers of the two components are the same, they 
have the same critical temperature and so may be ex- 
pected to begin condensation at the same time. However 
when the two components contain different number of 
atoms, the critical temperature of the larger component 
is higher {Tc ex N^^^), so it may be expected to begin 
condensation first. The time delay between the start of 
condensation of the two components will then depend on 
their atom number ratio and on the rate of evaporative 
cooling. We have performed simulations with large time 
delay such that one component begins to condense after 
the end of condensation of the second. We find that in 
this case the growth of the component which condenses 
later occurs at the edges of the first one, and they form a 
structure similar to the usual ground state, with no bub- 
ble pattern. This may be expected since the first compo- 
nent to condense creates a repulsive potential which re- 
pels the second component from condensing in the center 
of the trap. Thus, a sufficient temporal overlap between 
the condensation of the two components seems essential 
for formation of a bubble pattern. 



C. Comparison with experiment 

To make a concrete connection with the experimental 
results [1[, we simulated as close as possible the condi- 
tions under which the remarkable bubble patterns seen 
in Fig. (4b) and (4c) of Ref. [1[ were observed. The result 
is shown in Fig. Q . The number of atoms (reported in 
the caption) was also estimated from the experiment [18| , 
with 6 X 10"* ^^Rb atoms and 1 x lO" ^^Rb atoms. In this 
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FIG. 8: Dynamical growth of a dual BEG leading to the pattern formation seen in Fig. ((TJ. 



simulation we assume that the two components begin to 
condense together. Due to the different atom numbers 
this might not have been the case in the experiment, but 
should be good enough as long as the time delay between 
condensation of the two components is not too long. In 
the simulation We can see the formation of three bub- 
bles of *^Rb immersed in and between the larger cloud of 
^^Rb. This result is qualitatively similar to that observed 
in the experiment. The comparison with the experiment 
is quite encouraging, especially bearing in mind that we 
do not know parameters such as the exact time(s) of the 
growth of the two condensates, and the initial conden- 
sate "seed" sizes. The condensates in the experiment 
also show evidence for gravitational sag. In the radial 
direction, gravitational sag breaks the cylindrical sym- 
metry, an effect we have not attempted to simulate here. 
However we do include the effect of asymmetry in the 
axial direction due to residual gravity sag along this di- 
rection, which can occur due a slight tilting of the trap 
from the horizontal to the earth. 
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FIG. 9: The final spatial pattern resulting from the dynamical 
growth of a condensate with parameters approximating the 
experimental conditions of JJ], Fig. 4b there: ass = 81 bohr, 
ag7 = 99 bohr, ags-sr = 214 bohr. Trap frequencies are as in 
Fig. IT} above. The number of *^Rb atoms is 10* and that of 
^^'^Rb, 6 X lO**. The trap centers of the two components are 
shifted by 3.4 ^m from each other along the long axis of the 
trap, to simulate residual gravitational sag (see text). 
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IV. CONCLUSIONS 

In this work we have put forward a mechanism which 
explains the intriguing experimental observations of [l|. 
These observations showed the formation of a spatial 
modulated pattern of a dual BEC condensate upon the 
end of evaporative cooling in an highly elongated cigar 
trap. The mechanism we suggest includes two crucial el- 
ements. One, the growth of the dual condensates which is 
modeled by imaginary gain terms in the Gross-Pitaevskii 
gain equation, Eq. ([T]). Second, as the condensates grow 
and their mutual interaction becomes more and more 
dominant, the mechanism of modulation instability kicks 
in, and gives rise to the evolution of spatial modulation 
with smaller and smaller wavelength. The initial linear 
instabilities are amplified by the non-linearity in the sys- 
tem, and this process comes to an end when full spatial 
separation between the two components is achieved. The 
final pattern is very long lived due to the large energy 
barrier for the separated component "bubbles" to cross 
over each other. The quasi ID nature of the trap is es- 
sential, otherwise no barriers exist to prevent any formed 
"bubbles" from merging and coalescing into continuous, 
bulk components. 

Naturally, the gain model is a simplified, phenomeno- 
logical model of the growth of two coherent condensates 



from the thermal gases, which does not include the inter- 
actions between the condensates and the thermal gases, 
nor any direct temperature dependence. The thermal gas 
only enters through its role as a particle reservoir (for the 
gain terms). Thus, in principle, other mechanisms such 
as thermal fluctuations may also play a role in the final 
observed spatial pattern in the experiments. It would be 
of interest to have a more elaborate (but therefore, more 
complicated) modeling of the evaporative cooling process 
of a dual condensate. One such possible model may be 
the stochastic Gross-Pitaevskii equation [ij] which re- 
cently was shown to give excellent agreement with the 
experimental observation of spontaneous vortex forma- 
tion in evaporative cooling of a single component BEC 
|15J . Nevertheless, the encouraging agreement between 
our simulations and the experiment supports our belief 
that the simple model put forward here already captures 
the essential physics at the root of the observed phe- 
nomenon. 
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